Sum-rule Conserving Spectral Functions from the Numerical Renormalization Group 
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We show how spectral functions for quantum impurity models can be calculated very accurately 
using a complete set of "discarded" numerical renormalization group eigenstates, recently introduced 
by Anders and Schiller. The only approximation is to judiciously exploit energy scale separation. 
Our derivation avoids both the overcounting ambiguities and the single-shell approximation for the 
equilibrium density matrix prevalent in current methods, ensuring that relevant sum rules hold 
rigorously and spectral features at energies below the temperature can be described accurately. 
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Quantum impurity models describe a quantum system 
with a small number of discrete states, the "impurity", 
coupled to a continuous bath of fermionic or bosonic ex- 
citations. Such models are relevant for describing trans- 
port through quantum dots, for the treatment of corre- 
lated lattice models using dynamical mean field theory, 
or for the modelling of the dccohcrcncc of qubits. 

The impurity's dynamics in thermal equilibrium can be 
characterized by spectral functions of the type A BC (cu) = 
J ^e lLJt (B(t)C)T- Their Lchmann representation reads 



A BC (u) = J2{b\C\c 



a . b 



z 



(a\B\b) S(u - E ba ) , (1) 



with Z = e P Ea and Eb a = Eb~E a , which implies 
the sum rule J duiA BC (uj) = (BC)t- In this Letter, we 
describe a strategy for numerically calculating A BC (u) 
that, in contrast to previous methods, rigorously satisfies 
this sum rule and accurately describes both high and low 
frequencies, including uj < T, which we test by checking 
our results against exact Fermi liquid relations. 

Our work builds on Wilson's numerical renormaliza- 
tion group (NRG) method [1]. Wilson discretized the 
environmental spectrum on a logarithmic grid of ener- 
gies A - ", (with A > 1, 1 < n < N — > oo), with expo- 
nentially high resolution of low-energy excitations, and 
mapped the impurity model onto a "Wilson tight-binding 
chain" , with hopping matrix elements that decrease ex- 
ponentially as A - "/ 2 with site index n. Because of this 
separation of energy scales, the Hamiltonian can be di- 
agonalized iteratively: adding one site at a time, a new 
"shell" of eigenstates is constructed from the new site's 
states and the Mk lowest-lying eigenstates of the previ- 
ous shell (the so-called "kept" states), while "discarding" 
the rest. 

Subsequent authors [2, 3, 4, 5, 6, 7, 8, 9, 10] have 
shown that spectral functions such as A BC (uj) can be cal- 
culated via the Lehmann sum, using NRG states (kept 
and discarded) of those shells n for which ui ~ hr n / 2 . 
Though plausible on heuristic grounds, this strategy en- 
tails double-counting ambiguities [5] about how to com- 



bine data from successive shells. Patching schemes [9] 
for addressing such ambiguities involve arbitrariness. As 
a result, the relevant sum rule is not satisfied rigorously, 
with typical errors of a few percent. Also, the density 
matrix (DM) p = e~^ H jZ has hitherto been represented 
rather crudely using only the single iVr-th shell for which 
T ~ a.-H n t-i) [ g ] ; with a chain of length N = Nt ^ 

resulting in inaccurate spectral information for u) < T. 
In this Letter we avoid these problems by using in the 
Lehmann sum an approximate but complete set of eigen- 
states of the Hamiltonian, introduced recently by Anders 
and Schiller (AS) [11]. 

Wilson's truncation scheme. — The Wilson chain's ze- 
roth site represents the bare impurity Hamiltonian ho 
with a set of do impurity states \<7q). It is coupled to a 
fermionic chain, whose nth site (1 < n < N) represents 
a set of d states \cr n ), responsible for providing energy 
resolution to the spectrum at the scale A"™/ 2 . For a 
spinful fermionic band, for example, a n G {0, T, i, T-l}> 
hence d = 4. (Bosonic chains can be treated similarly 
[10].) The Hamiltonian H = Hn for the full chain is 
constructed iteratively by adding one site at a time, us- 
ing H n = H n —i + h n (acting in a cPdo-dimcnsional Fock 
space T n spanned by the basis states {\cr n ) ® ■ ■ ■ ® \a )}), 
where h n links sites n and n — 1 with hopping strength 

A~ n / 2 . Since the number of eigenstates of H n grows 
exponentially with n, Wilson proposed the following it- 
erative truncation scheme to numerically diagonalize the 
Hamiltonian: Let no be the last iteration for which a 
complete set {|s)^ } of "kept" eigenstates of H na can be 
calculated without trunction. For n > no, construct the 
orthonormal eigenstates {|s)^} of H n (the nth "shell"), 
with eigenvalues E™, as linear combinations of the kept 
eigenstates |s)^_i of H n -\ and the states \a n ) of site n, 



)n-l [Akx ] 



(2) 



with coefficients arranged into a matrix A£x whose el- 
ements are labelled by ss'. The superscript x = k or 
d indicates that the new shell has been partitioned into 



2 



matrix product states \s) n _^ 



,,["o+i] 



,|no+2) 



Kp+l) l°~no+2) K-l) 



discarded states \s') n 



-x 

^kk 1 ~ kept states \s') 



n 



FIG. 1: Diagram for the kept (or discarded) matrix product 
state \s')n (or |s')n) : the nth box represents the matrix block 
Akx ; its left, bottom and right legs carry the labels of the 
states |s)^_i, \a„) and |s')n (or |s')°), respectively. 



"kept" states (say the Mk lowest-lying cigenstatcs of H n ) 
to be retained for the next iteration, and "discarded" 
states (the remaining ones). Since h n acts as a weak 
perturbation (of relative size A" 1 / 2 ) on F n -i, the d-fold 
degeneracy of the states |er„) ® is lifted, resulting 

in a characteristic energy spacing A""/ 2 for shell n. It- 
erating until the spectrum of low-lying eigenvalues has 
reached a fixed point (for n = N, say), one generates 
a set of eigenstates {|s)^} with the structure of matrix 
product states [12] (Fig. 1). The states generated for the 
last, Nth shell will all be regarded as "discarded" [11]. 

Anders-Schiller basis. — Recently, AS have shown [11] 
that the discarded states can be used to build a complete 
basis for the whole Wilson chain: the states {|s)*} de- 
scribing the nth shell are supplemented by a set of d N ~ n 
degenerate "environmental" states {|e„) = \<tn) <£>••■ ® 
|<7„_|_i)} spanning the rest of the chain (|ej\r) denotes an 
"empty" state), to construct the set of states {|se)* = 
\e n ) ® |s)^}. These reside in the complete Fock space Tn 
of the full chain, spanning J-jq if n < no. Ignoring the 
degeneracy-lifting effect of the rest of the chain, these 
states become approximate eigenstates of the Hamilto- 
nian Hn of the full chain ("NRG-approximation"), 



H N \se)l ~ E™\se) 



(3) 



with eigenenergies independent of the (d N ~ n )-to\& degen- 
erate environmental index e„. (This will facilitate tracing 
out the environment below). By construction, we have 
°(se|s'e')° = 6 mn 5 e e > 6 SS > and 



0, 

5r 



KK 



..4 



KDj . 



m > n 
m < n. 



(4) 



The discarded states of shell n are orthogonal to the dis- 
carded states of any other shell, and to the kept states of 
that or any later shell. Combining the discarded states 
from all shells thus yields a complete set of NRG eigen- 
states of Hn, the "Anders-Schiller basis", that span the 
full Fock space Tn (%2 n henceforth stands for X^>n ) : 



= 5>>s, £M = EE i* e > 



(5) 



Local Operators. — Let us now consider a "local" op- 
erator B acting nontrivially only on sites up to uq. Two 



particularly useful representations are 



B = E Ml, WA, „ K ( S 'e| = e E 9£< • (6) 



The left equality, written B — B^k in brief, represents 
the operator in the complete basis set {|se)^ }, with ma- 
trix elements known exactly numerically (possibly up to 
fermionic minus signs depending on the environmental 
states; but these enter quadratically in correlation func- 
tions and hence cancel). The right-hand side (RHS) of 
Eq. (6) expresses B in the AS-basis and is obtained as 
follows: starting from Bk]°\ one iteratively refines the 
"kept-kept" part of B from, say, the (n— l)st iteration in 
terms of the NRG eigenstates {|se)*} of the next shell, 
including both kept and discarded states (x = k, d), 



Bt~ 1] = EEHn [4 n iL, r<*'ei = e*& 



(7) 



XX' ss'e 



thereby defining the operators S^J,, with matrix ele- 

ments [B^]^, = [A™ B& 1] A%)] Splitting off 
all xx' 7^ kk terms (dd, kd, dk) and iteratively refin- 
ing each kk term until n = N, we obtain the RHS 
of Eq. (6). It has two important features: First, the 
matrix elements of the time-dependent operator B(t) = 
e lHt Be~ lHt , evaluated within the NRG approximation, 
[B* x '(t)} ss , ~ [B™'] ss , e it( - E "~ E s'), contain differences of 
eigenenergies from the same shell only, i.e. calculated 
with the same level of accuracy. Second, by excluding 
kk terms it rigorously avoids the double-counting am- 
biguities and heuristic patching rules plaguing previous 
approaches [2, 3, 4, 5, 6, 7, 8, 9, 10]. 

Thermal averages. — To calculate thermal averages 
(••■)t = Tr[p...], we write the full density matrix 



P 



-0H 



jZ using the NRG approximation as in Eq. (3) 



EEH 



(se 



(8) 



where w n = d N - n ZfjZ and Z° = £° ^ m ■ The RHS 
of Eq. (8) expresses p as sum over , the density matrix 
for the discarded states of shell n, properly normalized as 
Tr [/Sdd] = 1, and entering with relative weight w n , with 
X)n w n = 1 • Similarly, for spectral functions we have 

(■•■)t = }]w n (...) n , A{u) = y]w n An(u;) , (9) 



where the averages (. . . )„ and spectral functions A n (uj) 

In] 

are calculated with respect to p D o of shell n only. 

Previous strategies [4, 5, G, 7, 8, 9, 10, 11] for ther- 
mal averaging amount to using a "single-shell approxima- 
tion" w n = & n N T for the density matrix and terminating 
the chain at a length N — Nt set by T ~ A-sWr-i), 



3 



As a result, spectral features on scales u> < T, which 
would require a longer chain, are described less accurately 
[see Figs. 2(a) and 2(b)]. Our novel approach avoids 
these problems by using the full density matrix (FDM), 
summed over all shells, letting the weighting function w n 
select the shells relevant for a given temperature yielding 
a smooth T dependence [see Fig. 2(c)]. Since w n has a 
peak width of five to ten shells depending on A, d and 
Mk and peaks at n-values somewhat above Nt [arrow 
Fig. 2(b)], spectral information from energies well below 
T is retained. 

Let us now consider the spectral function A bc (lo), for 
local operators B and C. Equations (4), (6), (8), and (9) 
can be used to evaluate (B(t)C) n . Fourier transforming 
the result we find (sums over ss' and a n implied): 



\m=n\ 



[m<nl 
Pkk 



n ^kk 

= E E 

m>n xx' 



\ r [m] [mn 
[4'x l J xx 



xx' 



= 5 SS * 



-0E" 



A Wm+l] 



S(co - EX) , 
(10) 



.[<x„] [nn] ,[<7„]t 4 [o- m +l]t 

^KD yUU ^IDK • • • %K 



Similarly, the static quantity (BC) n equals the first line's 
RHS without the S function. The matrix elements 
[AT 1 ]^, = E e m( se \P [ "M s ' e )m are given by the sec- 
ond and third lines, together with pfe = p£S <n ' = 0. 
After performing a "forward run" to generate all relevant 
NRG eigenenergies and matrix elements, A BC (u) can be 
calculated in a single "backward run" , performing a sum 
with the structure J2m>n ^P TCd ^ ' <"K)]'™'> starting from 

a j tt [ml. rod [mn] , , , , 

m = N. Here pW = Z^„> m w n p l xx (updated one 
site at a time during the backward run) is the full re- 
duced density matrix for shell m, obtained iteratively by 
tracing out all shells at smaller scales A - "/ 2 (n > m). 

Equations (8)-(10) are the main results of our "FDM- 
NRG" approach. They rigorously generalize Hofstetter's 
DM-NRG [8] (which leads to similar expressions, but 
using w n — S n N T and without excluding kk matrix el- 
ements), and provide a concise prescription, free from 
double counting ambiguities, for how to combine NRG 
data from different shells when calculating A bc (lu). The 
relevant sum rule is satisfied identically, since by con- 
struction J dujA BC (uj) = {BC) n holds for every n and 
arbitrary temperature and NRG-parameters A and Mk- 

Smoothing discrete data. — We obtain smooth curves 
for A BC (to) by broadening the discrete <5 functions in 
Eq. (10) using a broadening kernel that smoothly inter- 
polates from a log-Gaussian form (of width a) [2, 4] for 
\u>\ > ujq, to a regular Gaussian (of width u>q) for \ui\ < loq, 
where ujq is a "smearing parameter" whose significance is 
explained below. To obtain high quality data, we com- 
bine small choices of a with an average over N z slightly 
shifted discretizations [3] (see [13] for more details). 



Application to Anderson model. — We illustrate our 
method for the standard single-impurity Anderson model 
(SIAM). Its local Hamiltonian h = Eo- e o c o<T C Ocr + 



^ c OT c OT c oi c oi describes a localized state with energy eo, 
with a Coulomb penalty U for double occupancy. It 



is coupled to a Wilson chain J2 na \ n (cl l+la c nrT + h.c.), 
which generates a local level width T. We calcu- 
lated A<(uj) = A c °° c °°(-lu), A>(lu) = A c ^ c o,(lu) and 
A = A> + A < . An "improved" version A im thereof 
can be obtained by calculating the impurity self energy 
E(w,T) [6, 13] via FDM-NRG, which is less sensitive to 
smoothening details and yields more accurate results for 
the Kondo peak height *4t~o(0) at zero temperature. 

Sum rules. — As expected, we find FDM-NRG to be 
significantly more accurate at lower computational cost 
than NRG or DM-NRG [8, 15]. The sum rules 



dwA(u)) 



1 



(11) 



hold exactly to 10~ 15 for our discrete data, and to 10~ 4 
after smoothing (due to numerical integration inaccura- 
cies). Moreover, even for Mk as small as 256, our results 
for .4t~o(0) and A 1 ™^^) typically agree to within 2% 
and 0.2%, respectively, with the Friedel sum rule, which 
requires irTA^Q*" = sin 2 7r(co a c 0o .)o- The exact relation 
A < (uj) = f(ui)A(u>) (f is the Fermi function), which 
follows from detailed balance, is likewise satisfied well 
(though not rigorously so): the left-hand side of Eq. (11) 
typically equals / dwf(u>)A(u)) to better than 10 -4 . 

Low-frequency data. — Due to the underlying logarith- 
mic discretization, all NRG-bascd schemes for calculating 
finite-temperature spectral functions inevitably produce 
spurious oscillations at very low frequencies |cj| -C T. 
The scale St at which these set in can be understood as 
follows: the Lchmann sum in Eq. (1) is dominated by 
contributions from initial states \a) with energy E a ~ T, 
represented by NRG shells with n near Nt- The char- 
acteristic energy scale of these states limits the accuracy 
obtainable for energy differences E\, a to accessible final 
states 1 6). Thus the scale St is set by those shells which 
contribute with largest weight w n to the density matrix. 

We analyze this in more detail in Figs. 2(a) and 2(b) 
by purposefully choosing the smearing parameter to be 
unconventionally small, wo < T. The resulting spuri- 
ous oscillations are usually smeared out using luq > St 
[Fig. 2(a), thick gray (red) curve], resulting in quantita- 
tively accurate spectral functions only for \u\ > ujq ~ St- 
For conventional NRG approaches, the "single-shell" ap- 
proximation w n = S n N T typically leads to St — T, as 
can be seen in Fig. 2(a) [dashed (green) and thin solid 
(blue) line]. In contrast, FDM-NRG yields a signifi- 
cantly reduced value of St — T/5 [Fig. 2(a), black line, 
and Fig. 2(b)], since the weighting functions w n [inset of 
Fig. 2(b)] retain weight over several shells below Nt, so 
that lower-frequency information is included. 
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FIG. 2: (color online). FDM-NRG results for the spectral 
function At(u) of the SIAM, with U = 0.12, T = 0.01, e d = 
-U/2 (T K = 2.185 10~ 4 ), A = 1.7 and M K = 1024, unless 
indicated otherwise. Inset of (a): FDM-NRG result for At(u>) 
with uo in units of bandwith. For (a),(b) an unconventionally 
small smearing parameter was used, coo — 0.005 T (except 
for thick gray (red) curve in (a), with coq = 0.5 T), leading 
to spurious low-frequency oscillations. These illustrate the 
differences (a) between NRG (dashed green curve), DM-NRG 
(solid thin blue curve) and FDM-NRG (black curve) results 
for the regime to < T; and (b) between different choices of 
Mk and A for FDM-NRG, which yield different shapes for 
the weigths w n [shown in inset of (b)]: larger A reduces the 
scale St at which oscillations set in, but yields less accurate 
values for Kondo peak height in the regime St < to < Tk. 
(c),(d) Comparison of high quality FDM-NRG data (dots, 
solid curves) with exact Fermi-liquid results (black dashed 
lines) for (c) the conductance G(T) for T -C Tk; and (d) 
for At 1 (to) for T, to <C Tk. In (c), cg t was found from a 
data fit to c fit (r/r K ) 2 for T < T m (arrow). In (d) we plot 
5A T (to) = {AT (to) - A^iO^/A^iO) vs. lo/T k (curves) and 
SA(T) = [A^(0)/Aj, m (0)-1] vs. (T/T k )tt/VS (dots), for a set 
of 12 temperatures between 0.001 and 0.069 Tk (with curves 
and dots having same T in the same color), to illustrate the 
leading to- and T-behavior of At 1 {to); the dashed black line 
represents the expected Fermi liquid behavior in both cases, 
-TT%a; 2 vs. x. 

Fermi-liquid relations. — To illustrate the accuracy 
of our low-frequency results, we calculated A™ 1 ^) for 
u>,T <C Tk for the symmetric SIAM, and made quan- 
titative comparisons to the exact Fermi-liquid relations 

[141 



A T (cu) ~ A 



C ,T s2 3c , U ,2 

_ 2^Tk' ~ 2^Tk ' 



(12) 



duo A (oo, T) (- ~ A [l - c(|-) 2 ] . (13) 



Here Aq = l/7rr, c = 7r 4 /16, and the Kondo tempera- 



ture Tk is defined via the static magnetic susceptibility 
[4] Xo| T=0 = l/4T K . Figures 2(c) and 2(d) show the 
FDM-NRG data [gray (colored) dots and lines] to be in 
remarkable quantitative agreement with these relations 
[black dashed curves] . The results for the "conductance" 
G(T), being a frequency integrated quantity obtained by 
summing over discrete data directly without the need for 
broadening, are more accurate than for ^4™(w), and re- 
produce the prefactor c with an error consistently below 
5% (hitherto, errors of order 10-30% had been custom- 
ary). The smoothness of the data in Fig. 2(c), obtained 
using temperatures not confined to the logarithmic set 
[g ra y vertical lines in Fig. 2(c)], together with the 
remarkable stability with respect to different z-shifts il- 
lustrate the accuracy of our approach. 

Conclusions. — Our FDM-NRG method offers a trans- 
parent framework for the calculation of spectral functions 
of quantum impurity models, with much improved accu- 
racy at reduced computational cost. Its results satisfy 
frequency sum rules rigorously and give excellent agree- 
ment with other consistency checks such as the Friedel 
sum rule, detailed balance, or Fermi-liquid relations, in- 
cluding the regime uo < T . 

We thank F. Anders, R. Bulla, T. Costi, T. Hecht, W. 
Hofstctter A. Rosch and G. Zarand for discussions, and 
the KITP in Santa Barbara for its hospitality. The work 
was supported by DFG (SFB 631, De-730/3-1,3-2), and 
in part by the NSF (PHY99-07949). 

Note added. — Just before completion of this work we 
learned that Peters, Pruschke and Anders had followed 
up on the same idea [15]. 
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APPENDIX 

Smoothing Discrete Data:— The FDM-NRG yields the 
spectral function in the form of a Lehmann sum over dis- 
crete 5- functions, which have to be broadened to obtain 
a smooth function. Ideally, this should be done using a 
procedure for which the results are independent both of 
the parameters used for broadening and for discrctizing 
the conduction band. 

We calculate the smoothened spectral function us- 
ing A(uj) = j dbj'K (u,w')Aaw(w')- Here Aaw(w') = 
J2 n w n A n (uj') represents the discrete numerical data ob- 
tained from Eq. (10), which in practice we collect in 
binned form, typically using 250 bins per decade in uj' , 
so that J duj 1 becomes a sum over bins. The raw data 
is folded with the broadening kernel K(ui,u)'), which we 
choose to be of the following form: 

K(w,(j') = L(uj,uj')h(uj') + G(uj,uj')[l- h(uj')} (14a) 

where 



L(uj,uj = X- | v 



(14b) 



= ejuxS} ( igL^M +7 _ f ) 2 (7 _ f) 



vAFa|w '| e 



G(u>,ui') 
h{uj r ) 



{ 1 , 



(14c) 



i ^/^ V (14d) 
e v <* I , \ui | < uiq ■ 



The chosen kernel K constitutes a smooth interpola- 
tion, of somewhat arbitrary shape h(u>'), between a log- 
Gaussian [2, 4] broadening kernel L on the one hand, used 
for all ^'-frequencies but the smallest (with w and uj' re- 
stricted to have the same sign) ; and a Gaussian broaden- 
ing kernel G of width ujq on the other, used for \u>'\ < ujq 
to smoothly connect the regimes of positive and negative 
frequencies. We choose ujo to be roughly a factor of 2 
smaller than the smallest energy scale in the problem, 
including the Kondo temperature Tk (note that by con- 
struction in Eq. (14d) the transition to regular Gaussian 
sets in below ujq). 

The log-Gaussian kernel L(uj, uj') was purposefully cho- 
sen to have the following three desirable features: 

(i) Frequency- dependent width: being Gaussian on a log- 
arithmic scale, on a linear scale its width as function of uj 
is proportional to uj' . This is needed to deal with the fact 
that spectral data generated using Wilson's logarithmic 
discretization grid is more coarse-grained at large fre- 
quencies than at smaller ones. 

(ii) Conservation of weight: we have J du)L(u>,u)') = 1, 
ensuring that J duiA(u) = J duj'A Taw (uj'). 

(iii) Conservation of peak height: for the choice 7 = a/4 



(adopted henceforth) L is symmetric under uj <-> uj', so 
that also J doj'L(oj, uj') = 1. This ensures that the loga- 
rithmic broadening kernel maps a constant function onto 
itself (if A Ia , w (oj') = Aq, then A law (w) = Ao), and thus 
does not change the height of a peak whose width on a 
logarithmic scale is broader than a. (In this respect our 
L differs from that of [2, 4, 15].) 

Since choice (iii) implies that our log-Gaussian Kernel, 
as function of uj, describes a peak asymmetric w.r.t. uj' 
(shifted by a/4 on a log scale), on a linear w-scale the 
broadened data is stretched relative to the raw data by 
factor e a / 4 . This effect can be minimized by keeping a 
as small as possible. The smoothening of plain NRG data 
typically requires a ~ l/VA (e.g. 0.7 for A = 2). How- 
ever, smaller values (e.g. a < 0.3 or even smaller) can 
be achieved by using the "z-trick" [3] : collect several (say 
N z ) sets of discrete FDM-NRG data, each obtained from 
a different, slightly shifted logarithmic grid {A~™~ 2 } of 
discrete frequencies, for N z different values of z between 
—0.5 and 0.5, and average the results. The hopping ma- 
trix elements along the Wilson chain are recalculated for 
each z by carefully tridiagonalizing the underlying loga- 
rithmically discretized Hamiltonian. 

Self-energy representation: — The accuracy of the re- 
sults for A a [uS) for the Anderson model can be improved 
by expressing it in terms of the impurity self energy [G] : 
first, note that A a (uj) = —lm[Q^(uj)]/ir, where G^(uj) is 
the Fourier transform of G^'(t) = — i0(t)({c Oa (t), c' 0<7 })t- 
An improved version for Q^(uj) can be obtained by ex- 
pressing it as 
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^(uj) = U 



Here A a (uj) , the [/-independent part of the self energy 
which characterizes the level's broadening, can be com- 
puted exactly, G#(u>) is the standard ("non-improved") 
version of the correlator, and (uj) is the Fourier trans- 
form of -i6»(t)({[cJ ) _ CT c _ (T c 0(T ](t),c^) T . 

We calculate the imaginary parts of Q„(uj) and ^F^(uj) 
using FDM-NRG from Lehmann representations of the 
form (10), smoothen the discrete data as described above, 
Kramers-Kronig transform the smoothened results to ob- 
tain their real parts, and finally calculate E^(u>). Small 
wavy features in G&(uj) and ^^(uj) that reflect the loga- 
rithmic discretization grid largely cancel out in the ratio 
Yli (uj). Thus, smooth results for G™(uj) and A™(uj) can 
be obtained using much less (or even no) z-trick averag- 
ing, thus reducing the number of distinct FDM-NRG runs 
required to get good results. Moreover, since T%(uj — > 0) 
at T = is found to approach 0, the self-energy represen- 
tation also improves the accuracy with which the Fricdcl 
sum rule is fulfilled. 



